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Abstract: We investigate the thermal correlator of the trace of the energy- momentum 
tensor in the SU(3) Yang-Mills theory. Our goal is to constrain the spectral function 
in that channel, whose low-frequency part determines the bulk viscosity. We focus on the 
thermal modification of the spectral function, p{oj, T) — 0). Using the operator-product 
expansion we give the high-frequency behavior of this difference in terms of thermodynamic 
potentials. We take into account the presence of an exact delta function located at the 
origin, which had been missed in previous analyses. We then combine the bulk sum rule 
and a Monte-Carlo evaluation of the Euclidean correlator to determine the intervals of 
frequency where the spectral density is enhanced or depleted by thermal effects. We find 
evidence that the thermal spectral density is non-zero for frequencies below the scalar 
glueball mass m and is significantly depleted for m < w < 3m. 
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1. Introduction 

Understanding the emergence of bulk properties and collective behavior from the strong 
interaction of particle physics is a goal that has been pursued vigorously in recent years. 
Experimentally, heavy-ion collisions have provided a way to heat up nuclear matter to 
temperatures at which quarks and gluons are expected to dominate its macroscopic prop- 
erties [1, 2, 3, 4]. Theoretically, both analytic [5] and numerical methods [6] are being used 
to unravel the intricate dynamics of QCD at finite temperature and baryon density. 
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At temperatures of a few hundred MeV and small baryon density, lattice QCD remains 
the tool of choice to determine the equilibrium properties of the system. As we shall discuss 
in detail, access to the near-equilibrium properties, such as the transport properties, is 
limited in this approach, because lattice QCD employs the Euclidean formulation of thermal 
field theory. Real-time properties can thus only be determined by analytic continuation, 
which for numerical purposes represents an ill-posed problem. However, it is not excluded 
that we can learn about the gross features of the spectral density by suitably combining 
numerical and analytic methods. The purpose of this paper is to illustrate how such a 
combination of methods can be put to practice, taking as an example the bulk channel in 
the SU(3) Yang-Mills theory. 

One of the lasting legacies of the Relativistic Heavy Ion Collider (RHIC) experiments is 
the success of the hydrodynamic description of the collisions, see for instance the review [7]. 
The early agreement between ideal hydrodynamics and experiment was gradually refined to 
incorporate the dissipative effects of shear viscosity tj, the three-dimensional geometry [8, 9, 
10] and to estimate the sensitivity to initial conditions [11, 12]. Most recently, the effects of 
bulk viscosity ( have been investigated [13, 14]. At temperatures sufficiently high for weak 
coupling methods to be applicable, the bulk viscosity is three orders of magnitude smaller 
than the shear viscosity [15]. However, by analogy with many condensed matter systems 
(see for instance [16]), where a maxixum in the bulk viscosity is observed near liquid- 
gas phase transitions, one generally expects the bulk viscosity to become appreciable, if 
anywhere, near the rapid crossover between the confined and the deconfined phases. It 
has been suggested that this sudden onset of dissipation under expansion could trigger the 
freeze-out in heavy-ion collisions due to cavitation [14]; the breakup of a fluid occurs when 
a component of the stress-energy tensor, say T zz , turns negative. The behavior of £ near 
the QCD crossover was studied about two-and-a-half years ago, first using a QCD sum rule 
and a crude model for the bulk spectral function [17, 18], and subsequently using Euclidean 
correlators obtained by Monte-Carlo methods [19]. Both types of analyses missed the fact 
that a homogeneous fluctuation in T MM has a component which takes arbitrarily long to relax 
to equilibrium. It is therefore necessary to revisit these analyses, and this constitutes the 
second motivation for this work. A study based on AdS/CFT methods [20] also found an 
increase in the bulk viscosity near T c , however less pronounced than suggested by [17, 19]. 

In the near future, the low-energy runs at RHIC, and further ahead, the heavy-ion 
experiments at FAIR will explore the QCD diagram at larger net baryon density. It is 
widely thought, although not beyond doubt, that a chiral critical point exists in the T vs. 
Hb plane [21]. The associated critical exponents [22] indicate that the bulk viscosity should 
exhibit a strong divergence near the critical point. If this increase in £ is not too closely 
localized around T c , it could be one of the signatures of the critical point. 

The outline of this paper is as follows. We start in section 2 with a review of retarded 
and Euclidean correlators, establishing the precise relations between them. Details are 
collected in the appendix. In section (3) we turn to the specific case of the bulk channel 
in the SU(A^) gauge theory. We describe how to use the bulk sum rule as a constraint on 
the vacuum-subtracted spectral density Ap, and present the operator product expansion 
(OPE) prediction for the large-frequency behavior of the spectral density. In section 4 we 
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give a numerical application of our results for N = 3. We first present data deep in the 
confined phase, as well as in the deconfined phase. The sum rule is combined with lattice 
data to yield information on the changes in the spectral density that occur between the 
confined and the deconfined phases. We end with the conclusions of this study. 

2. A review of finite-temperature correlators 

In this section, we give the definitions of the retarded correlator and the Euclidean cor- 
relator, and show their relation in detail. This allows us to fix our conventions, but also 
to highlight a subtlety that arises in the zero-frequency case. The retarded correlator, 
evaluated at a purely imaginary frequency, is the analytic continuation of the Euclidean 
correlator, which is only available at the Matsubara frequencies. We show that the odd 
derivatives of the retarded correlator along the imaginary axis, evaluated at the origin, 
coincide with the derivatives of the spectral density at the origin. The even derivatives, 
on the other hand, correspond to integrals over all frequencies of the spectral density. We 
illustrate the relations that we derive on the example of the bulk channel in SU(iV) gauge 
theory. Working in frequency space rather than in time-coordinate space, we show the 
relation between the bulk viscosity and the frequency-space Euclidean correlator, Fig. (1). 
We work in the canonical ensemble, i.e. expectation values of operators are defined by 

(0) = I tr { e -^ O}, (2.1) 

where the partition function Z{f3) is such that (1) = 1 and f3 = 1/T is the inverse tempera- 
ture. In the following we consider a generic Hermitian operator O. We use the Heisenberg 
representation, with 

6{t) = e i6t 6(0)e~ iAt . (2.2) 
We start by introducing the real-time correlators, 

G R (t) = ie(t)([6(t),d(o)}}. (2.3) 

Its Fourier transform 

P oo 

G R (co) = i dte iut ([d(t),d(0)})- (2-4) 
J o 

plays an important role in the following. With these conventions we note that Gr is 
analytic in the half complex plane Im(u;) > 0. 
The Euclidean correlator is defined as 

G E (t) = (6(-it)O(0)). (2.5) 

It obeys the Kubo-Martin-Schwinger relation 

G E (P ~t) = G E {t) , (2.6) 
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and can be expressed as a Fourier series on the interval < t < (3, 

G E (t) =T^G i E ) e^ t , (2.7) 



eez 

Gf = / dte- iuJit G E {t), (2.8) 
J o 

where uj£ = com ■ £ = 2vrT £. 

The spectral representation of the correlators is useful to prove simple relations between 
the Euclidean and real-time correlators. One finds that (see the appendix) 

G R (iu, e ) = G$, i^O. (2.9) 

The relation of the Matsubara zero-mode £ = with the zero-frequency retarded correlator 
requires special care. We shall see that their difference can be expressed in terms of the 
spectral function. 

The spectral function is defined as 

p(w) = -ImG fl (w). (2.10) 

7T 

For real uj and in finite volume, p is a distribution, namely a discrete sum of delta functions. 
However, in the infinite-volume limit of an interacting theory, it is typically a smooth 
function, except possibly at a finite number of points. The most useful relations among 
the correlators are those that survive the infinite- volume limit. In the appendix, we show 
that 

<4 0) - G R {iui) = lim / ^p( w )^L^, Vu, 7 >0. (2.11) 

In the limit ojj — > 0, only the integration region around the origin contributes on the 
RHS. If the spectral function is finite at the origin, the RHS is 0(w/), and therefore 
lim e ^o Gji(ie) = G^\ On the other hand, if the spectral function contains a delta-function 
at the origin, p(uj)/uj = A6(u), lim e ^o Gji(ie) = G^ — A. We have written Eq. (2.11) with 
an explicit UV cutoff to make it clear that one should take the limit ujj — > first. 

The connection between the retarded and Euclidean correlators is thus precisely es- 
tablished in frequency space. In appendix we also rederive the relation between them in 
coordinate space (well-known in particular to lattice practitioners), 

[°° , x cosho;(f -t) 
G E (t) = l dup(w) sinh ^ /2 i . (2.12) 

which we shall use in section (4). 

2.1 The retarded correlator along the imaginary axis 

If we assume for a moment that 

roc j 

/ p(u) — , (C>0), (2.13) 
Jc w 
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converges in the UV, a simple calculation based on the explicit expression for the spectral 
density (A. 13), or alternatively a contour integration, then leads to the Kramers-Kronig 
relation 

dup(u) 2 , 2 = G R (iun) , Vwi > 0. (2.14) 

In general, it is necessary to perform a subtraction from p to make the dispersion integral 
converge in quantum field theory. Starting from (2.14), it is easy to derive by recursion 
the following equalities, 



j2n poo 

Mr^^/^W-V)^ (2.i5) 



(-i>" +1 -fSJre RM = rW MI W^i^./ w "(« = o). (2.16) 

Thus at the origin, the even derivatives of the retarded correlator along the imaginary axis 
are given by integrals over the spectral density, while the odd derivatives are equal to the 
corresponding derivatives of the spectral function. 

2.2 Remarks on the Euclidean correlator 

Combining Eq. (2.14) with Eq. (2.11), one obtains 



4°)=/°°^^. (2.17) 



Equations (2.11) and (2.17) show that while Gr{iuji) is insensitive to the presence of a 
uj8{u) term in the spectral function for any oj\ > 0, G^ definitely receives a contribution 
from it. 

Secondly we remark that the positivity of p{u)/u implies the positivity of certain linear 
combinations of the Euclidean correlator. For instance, 

4°>-<>^f/>^^, (, 18 > 

2.3 An example 

Here we anticipate the results of section (3) to illustrate the general formulas of this section. 
In Fi gure (1), we sketch a possible form of the retarded correlator for the operator 
projected onto zero spatial momentum. In order to deal with a UV finite function, we 
subtract the corresponding zero-temperature correlator, 

AG R {iou) = G R (iu, T) - G R (iio, T = 0) (2.20) 

and leave the temperature-dependence implicit. This subtraction is sufficient to make the 
dispersion integral (2.14) converge. As we shall see in section (3.1), the spectral function 
of the operator f ci 3 xT At/i (x) admits a delta function located at the origin, with a weight 
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Figure 1: A sketch of the retarded correlator for the operator L -3 / 2 J d 3 xT AlAl (a;) at purely imag- 
inary frequencies (dotted line). The large dots indicate the values of AG^) , calculable in the 
Euclidean theory. Figure (2) shows that AGn(ie) = 3(e+p)(l — 3c 2 ) — 4(e — 3p) is negative in the 
deconfined phase. At high frequencies, AGr is constrained by the OPE. 



given by a thermodynamic quantity. Applying Eq. (2.11), we obtain AGn(ie) = AG^ — 
Tc v (3c 2 — l) 2 . Furthermore an exact sum rule fixes AG^) in terms of thermodynamic 
functions, and AG^(ie) turns out to be negative. Using the operator-product expansion 
(section 3.2), one can show that AGjj is negative at large Matsubara frequencies, implying 
the same for AGji(iuj) at large co in view of Eq. (2.9). Finally, Eq. (2.16) and (3.3) show that 
{d / duj) AG n{iuj) at the origin is proportional to the bulk viscosity, which has to be positive 
by the second law of thermodynamics. All this information is summarized in Fig. (1). The 
behavior at intermediate frequencies oj is not known to us, in principle AGniioj) could 
even change sign in view of the subtraction made. It is also not known whether the plotted 
function ever exceeds unity for oj a multiple of 2nT. 

In a numerical Euclidean approach, the available information is Gn(ioJi), £ = 0, 1, . . . 
The low-frequency real-time properties of the system, such as the bulk viscosity and the 
associated relaxation time, are determined by the odd derivatives of Gr(iuj) at the origin, 
see Eq. (2.16). The latter are difficult to determine, because oj = is a non-analytic point 
of Gr. In this respect, Eq. (2.16) shows that if the spectral density has a gap m, as it 
does at zero temperature in a confining theory, then Gr(ioj) can be expanded in integer 
powers of oj 2 . Therefore subtracting the zero-temperature correlator only affects the even 
derivatives of AGn(ioj) at the origin. 

Working in frequency space has the advantage that the high-frequencies are decoupled 
from the discussion at this stage. There are, however, technical difficulties in determining 
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Euclidean correlators in frequency space. Unlike the t-dependent correlators, they are 
off-shell correlators, and therefore require additive counterterms in order to approach the 
continuum limit with corrections that are suppressed with a power of the lattice spacing. 
An example of this complication arises in the derivation of sum rules [23]. In this special 
case the sum rule precisely dictates what the counterterms are. In general, renormalization 
conditions have to be imposed in order to determine them. 

Finally we make the observation that a polynomial of degree n m oj going through the 
n + 1 points {Gn{ie)} U {G^}™ =1 , has a slope at the origin given by 

n! J. (-!)< G$-G R {ie) 

When expressed in terms of dispersion integrals as in Eq. (2.18-2.19), the degree of conver- 
gence in the UV of this expression improves with the order of the polynomial: with some 
combinatorics one can show that it is asymptotically p(co)/io n+3 for n even and p{u)/uj n+2 
for n odd. In particular, the slope at the origin of the n = 2 polynomial is given by 

-J-(3G R (ie) - 4G« + Gg>) = f° du P ^ 2 ) 2 , (2.22) 

2u M J-oo + e 2 (u 2 + u z M ) (Co 2 + Aw 2 M ) 

where the linear combination (2.19) appears. In the specific case of the bulk channel, the 
spectral integral (2.22) is UV-convergent without any subtraction to p(u). Thus the slope 
of the quadratic polynomial is negative-definite, just as the actual slope of the retarded 
correlator -^Gr{iuj). One might be tempted to use it as a Euclidean estimator for the bulk 
viscosity (. However, if the spectral function has a transport peak at the origin of width 
w o ^ wjw, expression (2.22) amounts parametrically to (wo/wjk, while ^Gr{iuj) ~ £. It is 
clear why that is: a polynomial interpolating between the Euclidean points is only a good 
approximation to the retarded correlator if the latter does not contain a scale ojq ojm 
over which it varies significantly. 

2.4 The reconstructed Euclidean correlator 

When analyzing Euclidean correlators obtained from Monte-Carlo simulations, it is in- 
structive to study in what way the Euclidean correlators differ from what they would be if 
the spectral density was unaffected by thermal effects [24]: 

Wms [ M) ^^l. (2.23, 
The identity (0 < t < p) 

sinh W /3/2 " ^ 6 [Z - M) 

' 1 m£Z 

allows us to derive in particular the exact relations 

G£ c (i,T;0) = J2G E (\t + mP\,T = Q), (2.25) 

meZ 

GT(t, T; \T) = G E (t, \T) + G E (f3 - t, \T) . (2.26) 
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These expressions are useful in practice when one is interested in a difference of spectral 
functions. We will exploit Eq. (2.25). Even if the zero-temperature Euclidean correlator can 
usually not be determined for arbitrarily large separations, in that regime it is completely 
dominated by the lightest state in the channel, so that its functional form is known and 
completely determined by one matrix element and an energy. The contribution to Gg c of 
a single state of energy E is easily obtained, 

6G E (t,T)~ ginhW2 . (2.27) 

Secondly, the correlator falls off very rapidly, so that the terms with large values of m in 
Eq. (2.25) make a very small contribution. 

3. The bulk channel in QCD 

In this section the goal is to provide the theoretical background for the analysis of the bulk 
channel in SU(iV) gauge theory. Thus we consider the operator 

Here /3(g) is the beta function, (3(g) = —bg 3 + . . ., b = and T^ v (x) is the energy- 

momentum tensor. Its spectral function pg(uj,T) and the difference 

Ape(oj,T)= Pe (oj,T)-p e (oj,0) (3.2) 

play an important role in the following. The spectral function determines the bulk viscosity 
through the Kubo formula 

( T = - lim ™±J. m 3.3 

In what follows we review the behavior that hydrodynamics predicts for the low- 
frequency part of the associated spectral function pg. It turns out that it contains a 
term ujS(uj). We then define an 'auxiliary' operator whose spectral function coincides with 
pg except for not having the delta function at the origin. Next we present the asymptotic 
large-w behavior of pg based on the operator-product expansion. Finally, the spectral 
function of the auxiliary operator obeys a sum rule given in section (3.2). 

3.1 Low- frequency prediction from hydrodynamics 

In this section we consider more general connected correlation functions of the energy- 
momentum tensor, 

G E „{t,K) = J ^ 3 xe ikx (r^(i,x)r pCT (0)) , (3.4) 

as well as their associated spectral function p^ up<J (u},]i,T) via (2.12). The spectral function 
for T33 (assuming the spatial momentum k is pointing in the direction 3) reads [25] 

P33,33(^,k,T) _e+p T s w 4 



w vr (oj 2 -c 2 s k 2 ) 2 + (T s u>k 2 ) 2 ' 



(3.5) 
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Now one shoud recall that for any smooth function such that f(x)dx = 1, |/(x/e) 
provides a representation of the delta function. Here we will exploit this fact for the 
function 

1 



Si?) 



ir(b 2 - 1) 



4 

X 



{x 2 — b 2 ) 2 + x 2 



1 



(3.6) 



and b = c s /(T s k) and e = T s k . In this way one finds that, in the sense of distributions, 
,. P33,33(w,k,T) e + p 2 

hm — ■ = r s + (e + p)ctd(oj) + ... (3.7) 

k->0 00 7T 

This is in contrast with what is obtained if one simply sets k = in Eq. (3.5), in which 
case one misses the delta function. In the shear channel on the other hand, 

Pl3 - 13(u; ' k ' T) = 1 ^ , 1 - ^S(oo) (3.8) 

and the delta function is suppressed by a power of k 2 . Thus combining the shear and sound 
channels, 

9 £J £J 3 W 7T 

Using the exact relations 

poo,oo(^, 0, T) = ^^wJCw) , (3.10) 

/9 00 ,fcfc(w,0,T) = -3(e + p)uoS(uo) , (3.11) 

one finds 

Ma^T) ^ p w , w (a,,0,T) = 9C e + p 2 _ 2 

£J £J 7T C 2 

When one constrains the spectral function using sum rules and Euclidean correlators, one 
is always dealing with integrals of the spectral function over all frequencies. Therefore it 
is useful to subtract the contribution of the delta function in (3.12). One way to achieve 
this is to work with the operator 

a = £-3/ 2 J d 3 x (T^(x) + (3c 2 - I)Tqo(x)) . (3.13) 

Its spectral function, which we denote by p*, satisfies 

p e (uo, 0, T) = p^u, 0, T) + ^(3c 2 - l) 2 co 6(u>) (3.14) 

c s 

and, in view of Eq. (3.12), is free of the delta function at the origin. At equilibrium, 
recalling that c 2 = ||, 

-\(0*)=P-& (3-15) 

and therefore its linear fluctuations are not correlated with those of the total energy 1 . This 
is the physical reason why a delta function at the origin is absent in its spectral function. 

1 Recall that we are working in the canonical ensemble. I thank G.D. Moore for pointing this out to me. 



-9- 



3.2 Large-frequency behavior from the operator product expansion 

In section (3.1) we have analyzed the low-frequency behavior of the bulk spectral function 
using the general hydro dynamic predictions. We now turn to an analysis of its asymptotic 
high-frequency behavior. The leading p# ~ a 2 s {uj)ui A behavior was given for instance in [26]. 
In the OPE language, this contribution corresponds to the insertion of the unit operator, 
and it is therefore temperature independent. Because it makes a large contribution to the 
coordinate-space Euclidean correlator, it is best to remove it by subtracting the vacuum 
correlator from the thermal correlator. We therefore need to know what the leading be- 
havior of the subtracted correlator and of the corresponding subtracted spectral function 
is. 

The operator product expansion of local gluonic operators was calculated a long time 
ago in the context of QCD sum rules [27], and the results of more recent calculations are 
summarized in [28]. In many cases, the Wilson coefficients are known at next-to-leading 
order. However, only the Lorentz scalar operators were kept on the right-hand side of the 
OPE, since the correlator was evaluated in the vacuum. Since we want to apply the OPE 
to a system at finite-temperature, the traceless part of the energy-momentum tensor, for 
instance, can also contribute. The OPE results therefore need to be generalized to the 
finite-temperature context [29] . We have reproduced the result of [29] for the bulk channel 
using the background field method combined with the Fock-Schwinger gauge [30], 

A^,0,T) ~ ^ ) (e-3p)--^ — J (e+p). (3.16) 

We remark that while the two-point function AG^(ia;) of the operator T^ u — \5^ v 9 con- 
tains a logarithmic divergence (which manifests itself in the Wilson coefficient of 9 being 
0(1/<7q ~ logl/o) in lattice regularization [23]), AGji(ioj) is UV-finite for the operator 9. 
Therefore, assuming the validity of the operator-product expansion, the only relevant scale 
in its Wilson coefficients is oj, and the dispersion relation (2.14) then leads to Eq. (3.16) [31]. 
In section (4), we will use this information to constrain the form of the spectral function 
at large frequencies. 

3.3 The Bulk Sum Rule 

The idea of constraining the form of finite-temperature spectral functions with sum rules is 
not new [32, 31, 17]. Here we focus on the bulk channel. In Euclidean space, the following 
sum rule holds [33], 

G%(T) - G^ e (T = 0) = T^ e -^ , (3.17) 

where e and p are the energy density and pressure of the finite-temperature system. 

We now convert identity (3.17) into a sum rule for the bulk spectral function with the 
goal to constrain the latter. Combining equations (3.17) and (2.17), we obtain 

r - Pe("M = T 5 d T e —^ . (3.18) 
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Next, we separate the contribution of the uj5(uj) term in pg from the rest of the spectral 
function. In other words, we want to turn Eq. (3.18) into a sum rule for p+, the spectral 
function of the operator (3.13), which is smooth at the origin. Using Eq. (3.14), one 
finds 

2 f°— Mw,r)-^( W ,0)] = 3(l-3^)(e + p)-4(e-3p). (3.19) 
Jo u 

A few remarks are in order. It would have been less straightforward to directly derive a sum 
rule for p+, because the Euclidean correlator of 0* contains a short-distance singularity. 
We remark that Eq. (3.18) was already obtained in [17], but the presence of the u8(oo) 
term in pg was missed. In [34], a finite spatial momentum k was used as an infrared 
regulator and the contribution of the zero-frequency delta function correctly identified for 
the first time. If one insists on expressing Eq. (3.19) in terms of pg, one should strictly 
write hm e _o JZo $T? A Pe(", T ) = 3(1 - 3c 2 s )(e + p) - 4(e - 3p). 

4. Numerical results from the lattice 

In this section we describe the lattice setup and the numerical results obtained by Monte- 
Carlo simulations. We employ the isotropic Wilson action [35], 

S S = \ll tr{l-P^(x)}, (4.1) 

where the 'plaquette' P^ u is the product of four link variables U^{x) £ SU(3) around an 
elementary cell in the (//, v) plane. As a discretization of L 3 ^ 2 Og, we use 

a 3^ e(x) = ^^^ Retr | Efc p ofc(x) + Efe< ^ [ p M(x) + p w(x + a 6)]} ) (4 .2) 

X X 

up to an irrelevant additive constant (its contribution vanishes in connected correlation 
functions). As a local update algorithm, we use the standard combination of heatbath 
and over-relaxation [36, 37, 38, 39] sweeps in a ratio increasing from 3 to 5 as the lattice 
spacing is decreased. The multi-level algorithm is used to reduce the variance of the 
correlator [40, 41]. To set the scale and to evaluate the lattice beta function dgQ 2 /dlog a we 
use the parametrization of the Sommer scale r$ ~ 0.5fm given in [42]. We denote the mass 
gap by m, given by the lightest scalar glueball, whose mass is 3.96(5)/ro = m = 5.30(8)T C . 
Finally we relate the notation used in the rest of this section to the general definitions of 
sections (2) and (3), 

Cg(t,T) = [3 5 G E ,g(t), C*(t,T) = /3 5 G £ ,,(t), (4.3) 
Cl ec (t,T;T Q ) = /3 5 G^(t,T;T ) , CT(i,T) = CT(t,T;0) , (4-4) 
Ap+{u,T) = p*(u,T) -p*(w,0). (4.5) 

The quantities with a * subscript are associated with the operator 0* defined in Eq. (3.13). 
The Euclidean correlators C... are dimensionless, renormalization-group invariant quanti- 
ties. The goal of this section is to constrain the function Ap + , the main results being 
Eq. (4.23) and Eq. (4.24). 
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Figure 2: The right-hand side of the bulk sum rule as a function of temperature. The spatial 
volume in the finite-temperature simulations is typically 28 3 . 



We start by describing a method we have used to evaluate the right-hand side of the 
bulk sum rule (3.19). We also evaluate the 'condensates' appearing in the OPE (3.16) of 
the bulk correlator. We obtain C* at finite temperature by subtracting the contribution 
of the 5 function. We then describe the properties of the correlator at temperatures deep 
in the confined phase. This correlator is subsequently used to produce a 'reconstructed' 
correlator based on (2.25). From then on we are directly dealing with integrals of the 
subtracted spectral function, A/?*. We reuse lattice data presented in [19], and combine it 
with the bulk sum rule to determine the gross features of the bulk spectral function. 

4.1 Right-hand side of the bulk sum rule 

When using the bulk sum rule (3.19), the question arises, how to best evaluate the combi- 
nation of thermodynamic quantities appearing on the right-hand side of the equation. Note 
that the enthalpy e + p and the conformality measure e — 3p are relatively straightforward 
to obtain from the expectation value of the electric and magnetic field strengths. What 
is then needed is a reliable way to obtain the velocity of sound c 2 s , or alternatively, the 
specific heat c v . 

The Euclidean sum rule (3.17) can be derived on the lattice. Here we will use that 
lattice sum rule [23] to determine the velocity of sound via 

dg Q z /dloga J e+p e + p ^ 

where the discretization of the trace anomaly operator is given by Eq. (4.2). We used 
the leading perturbative expression for the expression inside the large brackets, since based 
on that estimate it only departs from 4 by 0.008 or less. 
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The right-hand side of the bulk sum rule (3.19), as calculated on N T = 6 lattices, is 
displayed in units of e + p in Fig. (2) as a function of temperature. It is negative for a 
significant range of temperatures above T c , and its leading behavior in perturbation theory 
is [43] 

3(1 -^)-4l_| = -10.(^.(^) 3 + ... (4.7) 

It is therefore likely to be negative throughout the deconfined phase. We also note that for 
a s = 0.25, the right-hand side of (4.7) amounts to about -0.03, a very small value compared 
to those appearing in Fig. (2). 



4.2 Evaluation of the correlator C±{t,T) 

Since we want to work with the correlation functions of the operator O*, we have to evaluate 
its Euclidean correlation function C*{t). Lattice perturbation theory [44] and numerical 
studies [19] show that it is advantageous to work with the operator Og. To convert its 
correlation function to the correlation function C*, we use the equation (see 3.14 and 2.12) 

C e (t, T) = CS, T) + (3c 2 - I) 2 . (4.8) 

Here too a combination of thermodynamic quantities appears. We have obtained it by 
the same method as the right-hand side of the bulk sum rule, see the previous subsection. 
Close to the phase transition where c 2 is very small, the result is rather sensitive to the 
value of c 2 , since it appears in the denominator. That is the main reason we opted for 
a direct calculation of c 2 from the lattice sum rule (4.6), rather than parametrizing the 
pressure and obtaining c 2 by taking derivatives. In this way we avoid the dependence on a 
parametrization ansatz, which can be significant in a region where the derivatives are very 
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Figure 4: The 8 correlator at about half the deconfining temperature. The leading-order per- 
turbative prediction, with a s = 0.30, appears at small t. The contribution of the lightest scalar 
glueball is given at larger t. On the right panel, the correlator has been divided by the contribution 
of the lightest scalar glueball. Treelevel improvement has not been applied. 



large. On the other hand, at high temperatures 




and the difference between Cg and C* is parametrically negligible, since both are O(a^). 

The resulting correlators are displayed in Fig. (3). The lattice data for Cg is the same as 
in [19]. In particular, the lattice spacing is given by 1/oT = 8. At high temperatures, there 
is very little numerical difference between them, confirming the parametric expectation 
(4.9). On the other hand we observe that the difference between Cg and C* becomes 
appreciable very close to T c . At that point, C* is significantly smaller than Cg. 

4.3 Low-temperature correlators 

The Cg correlator at about T c /2 is displayed in Fig. (4). The agreement between the data 
obtained at two different lattice spacings suggests that the cutoff effects on the correlators 
are small. The dotted curve at small separation, which represents the leading perturbative 
prediction for the choice a s = 0.30, shows that by comparison the non-perturbative cor- 
relator is rather large. The right panel shows that starting at about 0.4fm, the lightest 
glueball represents a large contribution of the full correlator. The data is not accurate 
enough to tell how well that state saturates the correlator beyon 0.6fm. 

Although the finite size of the lattice implies that the temperature is about T c /2, 
rather than zero, it is very likely that the spectral function is already close to the zero- 
temperature spectral function, albeit in finite volume. This statement can be justified as 
follows. The confined phase can be thought of as a dilute gas of glueballs propagating 
freely, up to occasional collisions. Evidence in favor of this picture comes from a precision 
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thermodynamic calculation in the confined phase [45]. Thus, to leading order at low 
temperatures, 

— \ exp(-/3m) (4.10) 

c v « (pm)s(0) , (4.11) 

where m is the mass gap. Based on these approximate expressions, we find that the 
contribution of the exact delta function to the correlator Cq at t = (3/2 is 

[C e (t,T) - C(t,T)] t=m ~ (/3m) 2 (4.12) 

Numerically, at T = T c /2 this evaluates to about 0.006, while C e {T,(3/2) = 50(11). The 
area under the transport peak can be estimated using kinetic theory. The general expression 
is given in [25], 

£^-^>-<^f v " (4I3) 

where A is much larger than the inverse relaxation time, but lower than the thermal scale. 
At low temperatures this contribution is thus parametrically smaller than the contribution 
of the exact delta function, and it is indeed less than 15% of the latter at T c /2. 

In summary, there are strong reasons to believe that the bulk spectral function at 
T = T c /2 is close to the zero-temperature spectral function. Therefore it is legitimate 
to use these Euclidean correlators to obtain the reconstructed correlators at temperatures 
above T c , where they can be compared with the thermal correlators. 

4.4 Evaluation of the reconstructed correlator 

The 'reconstructed' Euclidean correlator is by definition what the Euclidean correlator 
would be if the spectral function remained unchanged with respect to the in vacuo spectral 
function. We want to evaluate this reconstructed Euclidean correlator (Eq. 2.23) in order 
to compare it with, and eventually subtract it from the thermal correlator. This evaluation 
is simplified by the fact that the 6 correlator falls off very fast at zero temperature. A first 
approximation is thus provided by (see Eq. 2.25) 

G T £ c (t, T; 0) « G E (t, T = 0) + G E ((3 - t, T = 0) . (4.14) 

The state of smallest energy contributing to the zero-temperature correlator is the lightest 
glueball. Its mass and the matrix element F$ = (vac|00|gluebail) have been calculated on 
the lattice (in this definition the norm of the states are set to unity). When /3m is large, 
we can therefore estimate the leading correction to Eq. (4.14), 

, FjcoshM(/V2-t) 

sinh/3M/2 ' V ' ; 

As examples we give the respective size of contributions (4.14) and (4.15) to the recon- 
structed correlator (in units of T 5 ) at a few temperatures of interest in the deconfined 
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phase, 



T = 1.02T C : 22.10(71) 0.081(11) (4.16) 
T = 1.24T C : 15.36(39) 0.117(16) (4.17) 
T = 1.65T C : 9.47(43) 0.132(23) . (4.18) 

As these examples show, (4.15) is a small correction to (4.14), in fact smaller than the 
latter 's statistical error. In subsequent figures we will therefore display (4.14) as a simple 
estimate of the reconstructed correlator. We note that expression (4.14) is a strict lower 
bound on the reconstructed correlator. 



4.5 Large-frequency contributions to spectral sums 

For a large value of co 2 , Eq. (3.16) implies to leading order, 



oc 



duj A . mx fllNa s {oj 2 )\ e-3p 3 fllNa s (u> 2 )\ 2 , 
Ap*(w,0,T) = - ( ) — ^ - 7 ( ) . (4.19) 



e + p J U2 uj \ Qtt ) e + p 4 \ 67r 

For instance, if we choose 0J2 = 12/ro ~ 4.8GeV, then a s (uu2) — 0.14 [46], and the second 
term on the RHS of Eq. (4.19) amounts to (-0.046). The first term then amounts to 

-0.11, - 0.053, - 0.025, - 0.011 (4.20) 

at 1.24, 1.65, 2.22 and 3.22T C respectively. This evaluation shows that the two terms on 
the RHS of Eq. (4.19) are about equally important at T = 2T C and that the second term 
becomes dominant for higher temperatures. 

The contribution of the tail of the spectral function to the Euclidean correlator at 
t = P/2 is given in leading approximation by 

J W2 T sinhu;/2T '* y ' V ' 

At T = 2T C , the exponential suppression is quite strong, e ~ w 2/2T ^ q q^s anc i Ap*(u>2) 
itself is of order (-0.02). This makes the large- frequency contribution to the Euclidean 
correlator negligible, and in any case much smaller than its current statistical error. 



4.6 The temperature dependence of the bulk spectral function 

Fig. (6) compares the thermal Euclidean correlators Cg(t,T) and C+(t,T) with the recon- 
structed correlator Cl ec (t,T;To). We see again that the difference between Cq and C* is 
small, except close to T c . However, the difference between the reconstructed correlator and 
is small too. Therefore the ^independent shift between C± and Cq is nonetheless crucial 
to correctly evaluate C±(t,T) — Cl ec (t,T). From Fig. (6), we learn in particular that the 
central value of this difference is positive at t = (3/2. While on a finer lattice one could 
reliably exploit more data points, here we will restrict ourselves to this middle-point, which 
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carries the most information on the low- frequency part of Ap*. Numerically we find, in 
terms of the subtracted spectral function (see Eq. (2.12) and Eq. (2.25)), 



X(T) = 



e + p 



X(T) 



1 



e + P Jo 



[C*(t,T)-Cl cc (t,T)] 



doj Ap + (w,r) 
T sinh 



t=p/2 > 

0.7 ±3.1 
0.26 ±0.13 
0.22 ±0.10 



On the other hand, the bulk sum rule (3.19) tells us that 



Y(T) = 



e + p 



o 



— Ap*(w,T) 



-0.928 ± 0.038 
-1.204 ±0.086 
-0.578 ± 0.099 



1.02T C 
1.24T C 
1.65T C . 



1.02T C 
1.24T C 
1.65T r 



(4.22) 



(4.23) 



(4.24) 



These two equations suggest an enhancement of the thermal spectral weight at low frequen- 
cies and a suppression at the higher frequencies. Indeed the higher frequencies contribute 
much more to the bulk sum rule than to Eq. (4.23). 

As a technical remark we note that since a large cancellation takes place when subtract- 
ing the reconstructed correlator, the quantity X(T) is susceptible to a large discretization 
error in addition to the quoted statistical error. In obtaining the results (4.23), we have used 
treelevel-improvement [44] separately for the finite-temperature and the zero-temperature 
correlators before calculating X(T). In [19] it was shown that after treelevel-improvement 
the discretization errors affecting Ce(t = (3/2, T) are mild. Also, we have ignored possible 
finite- volume effects, which could potentially induce a large relative error in X(T). New 
data on finer and physically larger lattices will be necessary to obtain X(T) with completely 
controlled systematic errors. 

4.7 Thermal rearrangement of the spectral weight 

We have collected a certain amount of information on the thermal spectral function rela- 
tive to its zero-temperature counterpart. Firstly, we know the asymptotic large-frequency 
behavior, thanks to the operator-product expansion. Secondly, two integrals over all fre- 
quencies of the subtracted spectral function are known. We also know that the difference 
A/?* is positive for u below the mass gap m, since the zero-temperature spectral function 
vanishes in that region, and we know the position and strength of the glueball peak. 

When determining properties of the spectral function, it is technically advantageous 
to separate the exact delta functions it contains from the contributions which are smooth 
in the infinite volume limit. Identifying the delta function at the origin in the bulk spectral 
function was an important step in this direction. Since we are interested in low-energy 
properties, it is also advantageous to work with a subtracted spectral function, such that 
the difference goes to zero at high frequencies. This is achieved by subtracting the zero- 
temperature spectral function. However, doing so one introduces, with negative weight, 
delta functions that correspond to the stable one-particle states at zero-temperature [47]. 
In the pure SU(iV) gauge theory, these are the glueballs. The position and weight of 
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Figure 5: A smooth parametrization of the thermal part of the subtracted spectral function, 
Ap*(w, T)/(u! ■ s), compatible with the bulk sum rule and the lattice correlator at t = f3/2. Recall 
that the ratio (/s is given by the intercept of this function. 



these delta functions can in principle be determined. They have been determined for the 
lightest glueball [48], but only rough estimates are known for the second (and last) stable 
scalar glueball. Therefore we will not attempt to systematically isolate the stable-glueball 
contributions in the spectral integrals X(T) and Y(T), but it is instructive to estimate the 
negative contribution of the lightest scalar glueball, 

_ 5 X(T) = F * s 1 = J L6 ( 2 ) L24T c U 25) 

G { ' T(e + p) sinhm/2T | 0.55(10) 1.65T C . V ' ; 



. ScY{ T) = -*5L_ = 3 '°W , ^ (4.26) 
V ; m(e+p) I 0.81(14) 1.65T C . v ' 



These contributions are quite large compared to X(T) and Y(T) themselves. We have seen 
that Y(T) is negative. Based on the OPE, only a small fraction of its magnitude can be 
explained by a depletion of the spectral density for a; > W2 = 12/ro ~ 4.8GeV. Therefore 
the bulk of the contribution to Y(T) must come from the region m < u < u)^- The depletion 
of the spectral density in the region of those glueballs that are stable or lie slightly above 
the two-particle threshold must therefore be the explanation for the sign and magnitude of 
Y{T). Whether the glueballs are completely 'melted' in the deconfined phase, in the sense 
that the spectral function shows no narrow peaks in this region, cannot be answered on 
the basis of the available information. 

While no detailed properties of Ap+ can be determined based on the limited information 
at hand, we want to show somewhat more quantitatively how Eq. (4.23) and (4.24) dictate 
the rearrangement of the spectral weight. For instance, in an interval of frequencies, one 
may ask whether f^ 2 ^j-Ap* is positive or negative. In order to address this question in a 
systematic way, we apply the following procedure (it is of course not unique). We select 
a binning in the frequency variable uj. We choose [0, wi] and [wi,^]- The upper end of 
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the second bin, 0J2, is chosen large enough for the OPE to be reliable. Since p*(a;,T = 0) 
vanishes for uj < m, it is of physical interest to choose uj\ on the order of m. 

Second, since we expect the (infinite- volume) spectral function to be smooth at least at 
high-frequencies, we parametrize Ap* by a second-order polynomial in uj 2 in each bin. We 
require continuity and differentiability at the boundaries between the bins. For uj > UJ2, we 
use the expression provided by the OPE. This then leaves 6 coefficients to be determined. 
Four equations are provided by the matching conditions at the edges of the bins, and two 
are provided by the known frequency integrals of Ap*, Eq. (4.23) and (4.24). Determining 
this simple parametrization of the spectral function thus amounts, in the present case, to 
inverting a 6 x 6 matrix. 

The result of this procedure is show in Fig. (5). The bands indicate the statistical 
error obtained by standard error propagation. These curves show qualitatively what was 
expected from the remarks below Eq. (4.24), namely the appearance of a positive spectral 
weight for uj < m, and a depletion relative to the vacuum spectral density for uj > m. 
Incidentally, this is qualitatively consistent with the prediction of perturbation theory, 
where the 0(a 2 ) low-frequency and O(a^) high-frequency contributions cancel eachother, 
resulting in a spectral sum / dujAp*/uj that is 0(a%) [34, 29]. 

We note that the values of the intercept in Fig. (5) at 1.24T C and 1.65T C , although 
strongly dependent on the fact that we have smoothened the spectral function, are in 
good agreement with our previous estimates for £/s [19], where the unsubtracted spectral 
function pg was directly parametrized. This may well be because in both cases we made a 
smooth ansatz for the spectral function, but in view of the large cancellation that occurs 
in Ap*, the consistency is not totally trivial. 

5. Conclusion 

We have collected information, obtained by analytic methods, on the spectral density 
of the bulk channel in SU(iV) gauge theory. We have reviewed the precise connection 
between the Euclidean and the retarded correlator in frequency space, illustrated in Fig. (1), 
recalculated the OPE of the bulk channel correlator, Eq. (3.16), and presented a new way, 
based on (Eq. 2.25), to compute the reconstructed correlator. The latter technique allows 
one to directly constrain the difference of thermal to vacuum spectral densities, Ap. 

Numerically, we have found that in the deconfined phase of SU(3) gauge theory, a 
significant depletion of the spectral density must take place in the interval m < uj < 
3m, while we have presented evidence that a non-vanishing spectral density appears at 
frequencies below the mass gap, obviously a finite-temperature effect. This qualitative 
conclusion comes from combining the bulk sum rule with numerical lattice correlators. 
With numerical data of higher accuracy, it would be possible to quantify more precisely 
how much spectral weight is clustered around the origin with frequencies uj <T. 

Unfortunately, the behavior of the bulk spectral function at low frequencies near 
the deconfining temperature remains largely undetermined. The Euclidean correlator 
C(t,T ~ T c ) clearly exhibits an enhancement compared to higher temperatures, even 
after subtraction of the i-independent contribution, see Eq. (4.8), which is quite large just 
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above the transition. However this enhancement can be explained to a large extent by the 
properties of the vacuum, as the comparison with the reconstructed correlator shows. At 
the moment we can state that a depletion of the spectral density has to take place for fre- 
quencies above the mass gap, in view of the bulk sum rule, but are unable at present to give 
a useful estimate of the spectral weight below the mass gap. This more careful reanalysis of 
the bulk channel has taught us that the effects of a potentially large bulk viscosity near T c 
are more subtle to detect in spectral integrals than previous work suggested [17, 19]. More 
accurate lattice data, in particular for the speed of sound near the phase transition, would 
allow one to determine the cumulated spectral weight below the mass gap. We believe that 
this goal can be achieved in the foreseeable future. 

Acknowledgments 

I thank K. Rajagopal and U. Wiedemann for interesting discussions on the bulk channel 
and G.D. Moore, D.T. Son, S.S. Gubser, S.S. Pufu and F.D. Rocha for a correspondence 
on the bulk spectral density. The simulations were carried out on a Blue Gene L rack and 
on desktops of the Laboratory for Nuclear Science, Massachusetts Institute of Technology. 
This work was supported in part by funds provided by the U.S. Department of Energy 
under cooperative research agreement DE-FG02-94ER40818. 

A. Spectral representation of finite-temperature correlators 

In this appendix we give the spectral representation of various correlators. We use this 
representation to obtain Eq. (2.11). 

A.l Real-time correlators 

Inserting two complete sets of energy eigenstates in the expression (2.3) for Gr{€), one 
obtains 




(A.l) 



m,n 



where we employ the notation 



E, 



"nm 



— E n — E,> 



O, 



(n\6(t = 0)\m) . 



(A.2) 



In the same representation, the Fourier transform (2.4) reads, for Imw > 0, 




(A.3) 



(A.4) 



In particular, for later reference, 




(A.5) 
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= ^\°nm\ 2 e-^ +E ^ 2 2 % 5i y7 /2) (^>0). (A.6) 
A. 2 Euclidean correlators 

The spectral representation of the Euclidean correlator (2.5) reads 

G E (t) = \Y. \Onm\ 2 e-^e E ^ , (A.7) 

n,m 

and after some algebraic manipulations, one also finds for the Fourier coefficients (2.8) 
r {i) _ P ST\ n i2 -p(E n +E m )/2 E lm sinh(/3£ ram /2) 

Z ^ ^ + ^nm P^nm/2 

Thus in particular we see that Eq. (2.9) holds. 

The relation of the Matsubara zero-mode i = with the zero-frequency retarded 
correlator requires special care. Using the spectral representation, one finds 

4 0) - = | E |Qn m | 2 e"^^^-)/ 2 Smh . ( ^7 2 /2) 2 "? . (A.9) 



In the next section, we shall see that the right-hand side can be expressed in terms of the 
spectral function. 

A. 3 The spectral function, and its relation to Euclidean correlators 

Taking the imaginary part of the retarded correlator, Eq. (A. 3), one obtains 

1mG R ( U ) = \Y. \Onm\ 2 e~^ ( 5 - ? ^ 2 ) (A.10) 

For a complex frequency just above the real axis, uj\ = e, we obtain 

lmG R {co + ie) = ^J2\°^rn\ 2 ^ En [S(u; + E nm )-6{cu-E nm )], oj G R. (A.ll) 

We have used one of the standard representations of the delta function, 

Symmetrizing the sum in Eq. (A.ll) with respect to m and n, one easily reaches the form 
(w G R) 

= - lim ImG fl ( W + i€) = -| sin ^ V |O nm | 2 e -^ +B ™)/ 2 <Hu, - £ nm ) . (A.13) 

We thus observe that the limit of Gr for co approaching the real axis is in general a 
distribution. 
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1. It is easy to verify, using Eq. (A. 7) and (A. 13), that the Euclidean correlator can be 
obtained from p via Eq. (2.12). 

2. Substituting (A. 13) into the right-hand side of Eq. (2.11), one recovers Eq. (A. 9), 
thus proving Eq. (2.11). 

3. A simple calculation based on expression (A. 13) formally leads to Eq. (2.14). 
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Thermal and Reconstructed Bulk Correlators at T=l .02T„ 
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Figure 6: Bulk-channel correlators at three temperatures. Here N T — 8 and treelevel-improvement 
has been applied. The open circles correspond to what the Euclidean correlator would be if the 
spectral density remained unchanged from a much lower temperature, To ~ T c /2. 
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